A novel general mapping for bosonic and fermionic operators in 

Fock space 

Alexej I. Streltsov*, Ofir E. Alon^, and Lorenz S. Cederbaunx'- 
Theoretical Chemistry, Institute of Physical Chemistry, Heidelberg University, 
Im Neuenheimer Feld 229, 69120 Heidelberg, Germany 
(Dated: October 14, 2009) 

Abstract 

In this paper we provide a novel and general way to construct the result of the action of any 
bosonic or fermionic operator represented in second quantized form on a state vector, without 
resorting to the matrix representation of operators and even to its elements. The new approach is 
based on our proposal to compactly enumerate the configurations (i.e., determinants for fermions, 
permanents for bosons) which are the elements of the state vector. This extremely simplifies the 
calculation of the action of an operator on a state vector. The computations of statical properties 
and of the evolution dynamics of a system become much more efficient and applications to systems 
made of more particles become feasible. Explicit formulations are given for spin-polarized fermionic 
systems and spinless bosonic systems, as well as to general (two-component) fermionic systems, 
two-component bosonic systems, and mixtures thereof. 
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I. INTRODUCTION 



The well-known time-dependent Schrodinger equation governs the dynamics of many- 
particle quantum systems in different fields of modern physics [1, 2, 3, 4, 5, 6, 7]. To define 
a quantum system means to specify its Hamiltonian. A quantum system is considered to be 
made of interacting constituent parts, usually treated as point particles with some known 
characteristics like mass and charge. If the interaction potential between these particles is 
known then the Hamiltonian of the quantum system is defined. Eigenvectors and eigenvalues 
of this Hamiltonian provide complete description of all the properties of the considered 
isolated quantum system. 

In general, exact solutions for many-particle Hamiltonians are not known and, therefore, 
different numerical approaches and techniques are in use . The basic 

and simplest approach is to represent the unknown many-body wavefunction as a linear 
combination of some known many-body wavefunctions, i.e., to expand the solution in a 
known basis set. To solve the problem means to find the expansion coefficients or their 
evolution depending on whether one performs time-independent or time-dependent studies. 

For a quantum system made of indistinguishable particles, e.g., fermions or bosons, along 
with the Hamiltonian one has also to specify the quantum statistics of the system. Statistics 
enters the many-body solution via the basis set used: if each and every basis function fulfills 
a defined quantum statistics then any linear combination of the basis functions also possesses 
the same statistics. From the other end statistics reduces the required size of the Hilbert 
subspace, i.e., it allows to operate with a smaller number of many-body basis functions. 

In this paper we specifically deal with a system of N identical particles. We utilize 
the commonly used many-body basis functions for bosons, permanents 5| constructed as a 
symmetrized product of N one-particle functions or orbitals, and for fermions, determinants 
which are antisymmetrized Hartree products of orbitals [8]. All many-body basis functions 
(determinants or permanents) are attributed to configurations where N particles reside in 
M orbitals; for bosons M > 1 for fermions M > N. For orthonormal many-body basis 
sets one can associate the number of the many-body basis functions used with the size of 
the respective Hilbert subspace. Intuitively it is clear that a larger number of independent 
many-body basis functions provides a better description of the many-body solution. The 
expansion spanned by all possible permutations for fixed M and N is referred to as a full 
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configurational expansion or a complete Fock subspace. If the many-body wavefunction is 
represented as a linear combination of permanents for bosons and determinants for fermions, 
then the Hamiltonian of the system as well as any other operator can also be expressed in 
second quantized form. 

Once the finite many-body basis set is specified one can construct the respective Hamil- 
tonian matrix. Straightforward implementation of the standard quantum mechanical rules 
requires construction and operations with this Hamiltonian matrix. For example, typical 
diagonalization or propagation schemes, like the short iterative Lanczos (SIL) [9], utilize 
as a standard building block the product of this matrix to a corresponding state vector, 
representing the many-body wavefunction. In this work we propose an absolutely different 
"ideology" which allows to get the required result of the action of the Hamiltonian on a state 
vector without construction of the Hamiltonian matrix at all. The proposed novel theory can 
be effectively and naturally applied to any general operator represented in second quantized 
form. The construction and operation with the corresponding matrix is not needed. We 
show the generality and applicability of this theory to systems of fermions, bosons as well 
as to the mixtures thereof. 

The proposed theory has already been implemented for bosonic systems 
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within 



the multiconfigurational time-dependent Hartree for bosons (MCTDHB) 12|. Applications 



to multi-boson long-time dynamics in double-well [10J] and triple-well llj traps have already 
been performed successfully. These studies demonstrate that the ideas proposed in this work 
are indeed very beneficial and can be considered as a breakthrough in computing dynamics. 

The structure of the paper is as follows. For the sake of expositional clarity we first 
consider in section [Til a system of spin-polarized fermions, i.e., the system where all spins 
of all fermions are identical, so each and every fermionic orbital can be occupied only by 
one particle. In subsection III Bl we enumerate and address the fermionic configurations in a 
simple, compact and practical way. The novel idea how to get the action of one- and two- 
body fermionic operators on a state vector as well as of the general Hamiltonian is explained 
in subsection III CI Subsection III Dl deals with other quantities like expectation values of 
the Hamiltonian and density matrices. In section II III we consider a system of bosons and 
apply the proposed theory to get the result of the action of a bosonic Hamiltonian on a state 
vector, without construction of the corresponding Hamiltonian matrix. In subsection IIII Bl 
we first describe how to map the fermionic onto the bosonic configurations. This allows 
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us to adopt, after small modification, the fermionic enumeration scheme for bosons. Next, 
in subsection IIII Cj we see how bosonic one- and two-body operators act on a state vector 
and then give explicitly the desired total action of the bosonic Hamiltonian on it as well as 
related quantities like expectation values of the Hamiltonian and other one- and two-body 
operators. In section [IV] we generalize our ideas and findings to multi-component systems. 
In subsection IIV Al we first deal with general fermions, i.e., the system where fermions 
with up and down spin projections coexist. Then we demonstrate the applicability of the 
theory to binary mixtures of bosons in subsection IIV Bl and to systems made of spinless 
bosons and spin-polarized fermions in subsection IIV CI Extension of the proposed theory to 
multi-component systems is given in subsection IIVDI In section |V] we discuss the practical 
implementation of the proposed ideas and section I VI I summarizes and concludes. 

II. THE CASE OF SPIN POLARIZED FERMIONS 

A. Hamiltonian of interacting systems and the state vector \SS?) 

Let us define the system first. We consider a general Hamiltonian in the second quanti- 
zation form with one-body and two-body interaction terms: 

H = J2 h kq b\b q + l -Y. W ksq ib{blhb q , (1) 

k,q k,s,q,l 

where the matrix elements hk q of the one- and Wk sq i of the two-body operators are assumed 
to be known. Three-body and higher-order interaction terms can also be included in an 
obvious way. To complete the definition of the system, one has to specify the commutation 
relations for the creation and annihilation operators b\ and b q . Here, we operate with the 
systems of indistinguishable fermions and, therefore, the usual anticommutation relations 
are fulfilled: b k b q + b\b k = S kq - 

We expand the generic state vector of the many-body system in a linear combination of 
Nconf known many-body basis functions \n): 

N t 

1 " conj 

|tf> = c * i") • ( 2 ) 

n 

Traditionally, Slater determinants are taken as |n) for fermionic systems. Using fermionic 
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creation operators b\ each Slater determinant is assembled as 

\n) = (biy 1 (bl) n2 ■ ■ ■ (b\ 4 y M \vac) , (3) 

where rij can be either "0" or "1" for spin-polarized fermions. n = (ni, 712, ris, • ■ ■ ,71m) 
represents the occupations of the orbitals that preserve the total number of particles n\ + 
n 2 + ^3 + ' ' ' + um = N, M is a number of the one-particle functions, here M > N, and 
I vac) is the vacuum. 

In the above written expansion Eq. ([2]) we did not specify explicitly the size N con f of the 
problem, i.e., the size or length of the vector. Let us now consider the configurational 
space spanned by all possible distributions of iV fermions over M fermionic orbitals, i.e., a 
complete Fock subspace of the respective configurational space. For spin-polarized fermions 
the size of such a complete Fock subspace 0, Q is 

AW=(^), (4) 

where (^) = fc! , ■ O n the other hand, N con f is the dimension, i.e., the number of the 
elements of any state vector of the system. 



B. Enumeration of the Slater determinants 



Now our goal is to provide simple and compact scheme for enumeration of the config- 
urations. More strictly, we can reformulate this enumeration as a requirement to map M 
integers nx, ^2, n^, ■ ■ • ,um characterizing each configuration, i.e., Slater determinant [see 
Eq. ([3])] to one integer addressing it as a coordinate (index) of the state vector. Clearly, 
there are many ways to solve this problem. Here we report on one that utilizes the so-called 
Combinadic numbers [jjj]. 

Every fermionic configuration can be represented as a vector with M components, filled 
by "1" or "0". The number of orbitals M must be larger than the number of spin polarized 
fermions N. Each appearing "1" means that the corresponding orbital is occupied by one 
fermion and "0" means that it is not occupied, i.e., it has a hole. Since the total number 
of "1" and "0" characterizing a configurational vector, i.e., its length, is M, there are N 
particles - "1" - and M v = M — N holes - "0". Here M v specifies the number of unoccupied, 
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i.e., virtual orbitals. The general fermionic configuration reads: 



■12 



'1 

I Tmmgimno - • • miiQ iiiiiii). (5) 

For example, for the system of N = 7 polarized fermions distributed over M = 10 fermionic 
orbitals there are N = 7 occupied orbitals and M v = M — N = 3 unoccupied ones. Then, for 
instance, the fermionic configuration 1 1011101011) means that the second, sixth and eighth 
fermionic orbitals are vacant, i.e., are not occupied, while the first, third, fourth, fifth, 
seventh, ninth and tenth orbitals are occupied by one fermion. The problem of enumeration 
of M dimensional vectors filled by "1" and "0" has been successfully solved in the context 
of Ref. [lij]; here we utilize it. 

For a given system with N fermions and M fermionic orbitals the address of every 
fermionic configuration can be uniquely defined either by specifying the positions of all 
the "l"s, i.e., particles, or alternatively by giving the positions of all holes, i.e., the "0"s. 
For the sake of definitiveness and without loss of generality we use the positions of the holes 
to specify the configurations. The general configuration defined in Eq. ([5]) is described by 
My = M — N holes placed at positions £1,1*2, • •• , *at„- ft is convenient to order the holes: 
%x < %i < £3 < • • • < im v . The address of this configuration is computed as follows: 

M v 

k=i x M v + l-k 

So, each and every fermionic configuration in a fermionic state vector |^) has its own unique 
address (index), defined by the numbers £1, £ 2 , • • • , im v - F° r the above considered fermionic 
configuration 1 1011101011) with N = 7 particles and M v = 3 holes, the "0"s are located at 
£1 = 2, £ 2 = 6, £3 = 8 when we count their positions from the left [see Eq. ([5])]. Using Eq. ([6]) 
we obtain the address of this configuration J(2, 6, 8) = ( 10 ~ 2 ) + ( 10 ~ 6 ) + ( 10 ~ 8 ) + 1 = 65. We 



if ■ ■ \ 1 1 (N + M v — i k \ 



note that the enumeration scheme provided in Ref. [14j uses positions of the "l"s counted 
from the right, while in our scheme we count the positions of the "0"s, i.e., holes from the 
left. The inverse problem, i.e., a restoration of the "holes positions" £1, £ 2 , ■ ■ • , im v according 
to a given address J also can be solved [l^] if needed. Our goal is fulfilled - spin-polarized 
fermionic configurations are enumerated in an easy and compact way. 

Let us summarize, N spin-polarized fermions distributed over M fermionic orbitals span 
the complete subspace of the Fock space of N conf = ( N+ ™ v ) configurations, here M v = M-N 



is the number of unoccupied fermionic orbitals. The dimension of any state vector of 
this system is N con f. Every fermionic configuration Eq. (jSJ) in the respective Fock subspace 
is characterized by the positions of M v holes placed at ii,^, ••• ,Hi v - We can attribute 
a unique address J to each fermionic configuration according to the rule Eq. J = 
Jih,^, • • • i^m v )- This Combinadic-based mapping is a one-to-one mapping. In particular, 
the index J takes all values between 1 and N con f. 

Now a general |\&), Eq. (j2J), can be rewritten in a specific form: 



where it is explicitly stated that every fermionic configuration, i.e., Slater determinant 
|J(i)) is specified by M v holes placed at (ii, *2j ■ • ■ ^m v ) = i- It has a unique address 
J = J(ii,i2, ■ ■ ■ , im v ) which can be computed according to Eq. (jSJ). The index (address) J 
runs over all N con j configurations. 

C. Applying operators to 

The Hamiltonian ([1]) is defined as a sum of terms, each of them is a product of creation 
and annihilation operators - a pair b\b q for the one-body and a quartet b\b\bib q for the two- 
body terms, scaled by respective prefactors, i.e., integrals h kq and Wk sq i- In principle, any 
other operator can be represented in a similar way as a sum of contribution of one-, two-, or 
higher-oder terms. The idea is quite simple: if one would know the result of action of each 
of these terms on a state vector \$), the sum of all of them would give the required total 
result of the action of the Hamiltonian on a state vector. Here we would like to recall that 
the total number of terms in the Hamiltonian is defined by the number of the orbitals used. 
In the most general case of M orbitals the total number of the one-body terms is M 2 , and 
the total number of the two-body terms is M 4 . Any symmetry in the problem, including 
hermicity, reduces these numbers. 

1. The action of one- and two-body operators 

Let us consider the action of a pair b\b q of creation and annihilation operators on every 
fermionic configuration Eq. (jSj) of a general state vector \^f). As a first step we consider a 
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j=i 



7 



specific, say b\b 3 term, that kills a particle in the third orbital and creates a particle in the 
second one. For spin-polarized fermions due to Fermi-Dirac statistics the term b\b 3 provides 
a non-zero action only on the subset of configurations |nj.,0, 1, • • • , Um) having n 2 = 0, 
n 3 = 1: 

&J,6 3 |ni,0, 1, • • • ,n M ) = (-l)°|m> !>0, • • • ,n M ), 
where can be "0" or "1". Similarly for a general b\b q case, n k = 0, n q = 1: 

6t 6?K n 2 , „ 3 , • • ■ , n«> = (-l)^(^) - K n 2 , • • • , n q - 1, • • • , n k + 1, - • • , n M >, 

where (— l)^-^*.?) ni is a prefactor ensuring correct fermionic statistics of the antisymmetrized 
wavefunction and the summation i e (k, q) runs over all occupations rii between the k-th and 
g-th orbitals. We can interpret this well known result as follows: operation of any even combi- 
nation of creation and annihilation operators on a configuration (determinant) results in the 
re-addressing of this configuration (determinant) to another one. Since the occupation num- 
bers of the incoming |n 1; n 2 , n 3 , ■ ■ ■ , tim) and resulting \n\, n 2 , • • • , n q — 1, • • • , n k +l, • • • , tim) 
configurations are explicitly known, the numbers of the orbitals with zero occupations, i.e., 
the holes positions are also available and therefore, according to Eq. (Q, we can also com- 
pute their addresses in a state vector Eq. ([7]). Applying the enumeration scheme (mapping) 
introduced above one gets: 

blb q \J(i)} = (-l)^\m},kei, q ^l 

In other words, the J-th configuration is re-addressed to a new configuration with index 
J' = J(i') with some sign prefactor. For a given configuration |J(i)) and fixed k, q the 
integer d k f is equal to the number of fermions, i.e., "l"s located between the k-th and g-th 
orbitals: d k f = Yliefk q) n * - Note that there is no difference in whether to count the number of 
fermions between the /c-th and g-th orbitals in the original configuration J or in the resulting 
J' configuration, namely, d k f = d k j,. It is important to stress that the operation with even 
combinations of creation and annihilation operators results in a single-valued re-addressing, 
i.e., an initial configuration having address J is re-addressed to a single configuration with 
address J'. 

Having at hand the result of the action of the pair b\b q on a general configuration (de- 
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terminant) we can find out their action on the total state vector [see Eq. (JTj)]: 




**) = b\b q |*> = Cjb\b q | J) = y, c j 9 w)) > 



(8) 



j=i j=i 




where at the last step we have changed variables of the summation index from J' to J. Note, 
that this closed-form result of the action of the basic one-body operator on a state vector 
has been obtained without referring to the matrix representation of the respective operator. 



J is obtained as a product of the Cjk q element of the incoming vector having address J kq 
scaled by the (— l) d j q fermionic prefactor. The configuration J kq is related to J by making 
a hole at k-th and filling a hole at q-th orbitals. The b k b q term acts only if the index k 
coincides with one of the holes' positions (ii, ■ ■ ■ ,k,--- , %m v ) and the index q with neither 
of them. If these two conditions are not fulfilled the respective configuration does not 
contribute at all to the resulting state vector. The addresses J kq = J(ii, ■ ■ ■ ,k,--- ,imJ 
and J = J {ii-, ■ ■ ■ ,<?,••• > i>M v ) are computed by Eq. (jSj). In Eq. (jSJ) the summation runs over 
the index J implying that for a given J one has to restore the holes' positions i\, • • ■ , iu v 
first. However, this complication can be easily avoided if, instead, one starts M v nested loops 
running over the positions of the holes i\ < i% < ■ ■ ■ < %m v directly. In the latter case all 
ii, 4 ■ • , ihi v are explicitly available, as well as the resulting address J = J(ii, ■ ■ ■ ,<?,•'' > UfJ 
via Eq. (jSj). Hence, at each step to get the desired element Cj q of the resulting vector 
one has to compute one integer dj q and apply Eq. flS]) only once to get the index J kq = 
J {iii ' ' ' ) k, • ■ ■ ,im v ) of the respective incoming configuration and its value Cjk q . 

The action of a general (k 7^ s 7^ I 7^ q) two-body b k b\bib q term on the incoming state 
vector can be obtained using a similar strategy: 

N t 

1 ' conj 

\* kslq ) = blb%b q I vl/) = b% I**) = Y C J Sk I . ( 9 ) 

j=i 



How to understand this result? Every element Cj q of the resulting vector having address 




To get the J-th element C kslq of the resulting state vector one has to take the element 
Cjksiq of the incoming vector having address J kslq and multiply it by (— l) d J sl+dSj prefactor. 
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Using Eq. (jSJ) we compute the addresses of the resulting J = J(ii, ■ ■ ■ , I, • ■ ■ , «mJ 

and incoming J ksl i = J(i 1; • ■ • , k, ■ ■ ■ , s, ■ • ■ , im„) configurations. The configuration J ksql is 
related to J by making holes at fc-th and s-th and filling holes at q-th and Z-th orbitals. The 
b' k bg term acts only if the index k coincides with one of the holes' positions (ij, ■ • • ,k,--- , %m v ) 
and the index q with neither of them. For a given quartet k, s, I, q every component J of 
the incoming state vector is characterized by the integer d k j sl + d s j computed as a sum of the 
number d s j of fermions located between the s-th and Z-th orbitals of the configuration J and 
d k j al - between the fc-th and q-th ones of the configuration J sl . The b\b\bib q term provides 
non-zero action only on the smaller subset of the configurations having n& = n s = 0, 
rii = n q = 1. Therefore, in practical computations one needs to address only this subset. 



2. The action of the Hamiltonian 

Considering configurations as coordinates of the state vector \^f), we have seen that the 
action of each term of the Hamiltonian on the state vector re-addresses the coordinates 
of the original state vector, multiplying them by some known prefactors. So, instead of 
constructing the full Hamiltonian matrix and performing matrix-to-vector multiplications, 
we obtained the same result by reordering the components of the incoming state vector 
according to the action of every b\b q and b\b\bib q term and by multiplying them by the 
corresponding integrals h kq and Wk sq i and summing the results up. 

The general Hamiltonian Eq. is a sum of the one-body h = J2k q ^kqKPq anc ^ two-body 
W = \ J2k s q i Wksqib\.blb[b q terms. Using the results of the previous subsection we find the 
total action of all one-body terms on the initial state vector 



k,q k,q J=l 

Cj = h kq Cj q . 

k,q 

Here, all the components Cj Q have been computed above using Eq. (jSJ). Analogously, we 
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sum up all the contributions from all two-body terms: 

N i 

k,s,q,l k,s,q,l J=l 

k,s,q,l 

Finally, the desired action of the Hamiltonian on the state vector is a sum of the two 
obtained above vectors in Eqs. fllOljl 1H : 

N i 

H\i&) = h\V) + W\V) = £ c j l J (i)>> ( 12 ) 

j=i 

Cj = ^ + cf, 



which fulfills our initial goal. This closed-form result of the action of the Hamiltonian on 
a state vector is constructed without building the respective Hamiltonian matrix or even 
without referring to its matrix elements. 



blbtkb q |*> =-£ W ksql |***> = e cf W)) . (11) 



D. Other quantities of interest 

In the preceding subsections we have successfully derived a simple and straightforward 
technique to operate with state vectors without representing the respective operators in 
matrix form. In particular, we have utilized the basic fact that the result of action of any 
operator represented as a sum of products of creation and annihilation operators on a state 
vector is equal to the sum of action of each of these terms. So, the action of any operator on 
a state vector (called incoming state vector) results in a new state vector (called resulting 
state vector). Let us now show that the expectation value of the respective operator can be 
immediately computed as a dot-product of the incoming and resulting state vectors. 



Indeed, according to the standard definition 15l I la. Il7|. for a given state vector |\&) the 



elements of the reduced one-body density matrix read: 

(¥| b{b g |¥> = (*| \b{b q = (*|**«) = £ C}C^ 



N, 



Pkq 



conj' 



(13) 



j=i 



where we substitute the result for b\b q from Eq. ([8]). Thus, the elements of the reduced 
one-body density matrix p kq are obtained as a dot-product of the incoming |^) and resulting 
l^ 9 ) state vectors. 
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Similarly, the elements of the reduced two-body density matrix are obtained as a dot- 
product of the incoming |^) and resulting [see Eq. (jUJ)] \^ kslq ^ state vectors: 

J * con/ 

Pksi q = <*| 6l6l6i^ |*> = <*| \b\b%b q m] = <*|*^> = £ C^C** (14) 



j=i 



So, in the discussed scheme the elements of the reduced one- and two-body density matrices 
appear very naturally. Moreover, the elements of the reduced three- and higher- order 
density matrices can be obtained in a very similar way. 

Finally, taking the result of the action of the Hamiltonian on an initial state vector 
H\ty) from Eq. (|12|) . we compute the respective expectation value of the Hamiltonian as a 
dot-product as well: 



N, 



(¥|#|¥) = <tt| H\*}\ = Y,C}Cf. 



conf 



(15) 



We have prescribe here a few expectation values of particular interest in a many-body theory. 
Other quantities can be represented and computed in a similar way. 



III. THE CASE OF STUCTURELESS BOSONS 



A. General remarks 



Here, we deal with a system of N identical interacting bosons. The terms structureless 
or spinless are often used to specify that the bosons do not have an internal structure. 
Our goal is to show that our ideas how to operate with state vectors without resorting 
to the matrix representation of the respective operators proposed above for fermions can 
naturally be extended and applied to bosonic systems. The generic Hamiltonian Eq. (pQ) 
introduced above for the system of spin-polarized fermions is also applicable to the system 
of structureless bosons with the only modification concerning the commutation relations of 
creation and annihilation operators. Here, we operate with the systems of indistinguishable 
bosons and, therefore, we use the usual commutation relations: bkb^ — b^bk = Sk q - 

The generic bosonic state vector is expanded as a linear combination of N con f known 
many-body basis functions \n): 

N t 

1 * conj 

|*)= ^Ca|fO, (16) 

n 
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where \n) are permanents which are assembled as: 



n) 



1 



(b\) ni (biy ■ ■ ■ (b^y M \vac) 



(17) 



v / ni!n 2 !n 3 ! • • -n M ] - 



n = (rii,n2,n3, ■ ■ ■ ,hm) represents the occupations of the orbitals that preserve the total 
number of particles rti + n 2 + n 3 + ■ ■ • + tim = N, \vac) is the vacuum and M is the number 
of the one-particle functions; for bosons M > 1. 

The number of elements N con f in the bosonic state vector is equal to the number of 
the configurations used. Here we consider the configurational space spanned by all possible 
distributions of N bosons over M orbitals, i.e., a complete subspace of the respective config- 
urational space. We recall that one of the key consequences of this completeness is that an 
action of any operator on the state vector results in a new state vector defined in the same 
configurational subspace. 

n 

The size of this complete Fock subspace is [13J: 



This number is equal to the size of the configurational space spanned by N spin-polarized 
fermions distributed over M' = N + M — 1 fermionic orbitals, see section III AL Therefore, 
there exists one-to-one mapping between the configurational spaces of the iV-boson system 
distributed over M bosonic orbitals and a fermionic system made of N fermions distributed 
over M' = N + M — 1 fermionic orbitals and vise versa. In other words, every fermionic 
configuration having N + M — 1 components (occupation numbers) should be attributed to 
a bosonic configuration characterized by M components. Let us compare these isomorphic 
bosonic and fermionic systems. What is also equal in these two systems apart from the total 
number N of particles? The number of the fermionic unoccupied orbitals M v = M' — N, 
i.e., the number of fermionic holes is equal to the maximal number of the bosonic holes: 
M v = M — 1. We note that the maximal number of bosonic holes appear in configurations 
like \N, 0, • • • ,0), |0, N, 0, • • • ,0), etc. Since we already know how to enumerate fermionic 
configurations in terms of holes, see section IHBl we can adopt the fermionic enumeration 
scheme to the respective bosonic system. 




(18) 
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Enumeration index J 


N Fermions, N + M — 1 orbit als 


N Bosons, M orbitals 


1 


111 nnn nr\\ 

| 1111 • • • 111 (JUL) • • • UU) 


\N, 0, ••• ,0) 




W M-l 


M-l 


2 


1 1111- •• 11 01 000 •••00) 






iV-1 M-2 


M—2 


J 


1 111110 i_n,o lino • • • o nioo nmi) 

5 3 4 3 6 


|5,3,4, - - - ,3,0,6) 


Nconf — 1 


| 000 - • 00 10 1111 • • • 111) 


|0 i _^0 J l,iV-l) 




M-2 N-1 


M-2 


Nconf 


| POO • • • 00 1111 • • • 111) 


1 0, - - - ,0,N) 




M-l N 


M-l 



TABLE I: Mapping between fermionic and bosonic configurations made of iV particles. The 
number of the unoccupied fermionic orbitals is equal to the total number of the bosonic orbitals 
minus one: M v = M — 1. Now, we can adopt the fermionic enumeration scheme for bosons. 

B. Mapping and enumeration of the permanents 

Our goal is to provide simple and compact scheme for enumeration of the bosonic con- 
figurations. To utilize the formal isomorphism between the bosonic and fermionic systems 
considered above we first have to show how to map, or attribute a fermionic configuration to 
the bosonic one. The rule is very simple: the number of bosons residing in the first bosonic 
orbital n\ is equal to the number of fermions occupying successively the lowest fermionic 
orbitals from the bottom till the first fermionic hole. In other words, the occupation of the 
first bosonic orbital is denned as the number of "l"s appearing in the (N+M— l)-component 
fermionic vector [Eq. §5§] up to the first hole, i.e., the first "0", when counting from the left. 
The occupation number of the second bosonic orbital n<i is defined as the number of "l"s 
between the first and second "0", i.e., between the first and second fermionic holes. The 
third bosonic occupation number is defined as the number of "l"s between the second 
and third "0" , and so on. Table fl] illustrates this mapping. 

Consequently, the bosonic occupation numbers and the positions of the fermionic holes 
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are simply connected: 

ni = h- 1, 

n 2 = i 2 -h-l, 

n k = (19) 

TlM-l = Hl-l — Hl-2 — 1, 

n M = N + M -i M -i - 1. 
So, having a set of the bosonic occupation numbers |ni, n 2 , ■ • • , %/) using this scheme we 
can restore the positions of the M — 1 fermionic holes (ii,i 2 , • • • , «m-i) and the respective 
fermionic configuration and vise versa. Now we can explicitly use the relations Eq. ffl9|) 
between bosonic occupation numbers and fermionic holes to uniquely address bosonic con- 
figurations utilizing Eq. (JOJ. Finally, omitting all intermediate steps we get the address of a 
generic bosonic configuration \ni,n 2 ,n 3 , ■ ■ ■ ,n M ) in the bosonic state vector: 

^ (N + M - 1 - k - T* , riA , . 

J(ni,n 2 ,...,n il f_i,n M ) = l+X) ( M — k )' (20) 

fc=i ^ ^ 

Here nM enters the expression implicitly via the identity N = m + n 2 + ^3 + ■ • • + tim- 

Let us summarize. The system of N bosons and M bosonic orbitals spans the complete 
subspace of the Fock space of N con f configurations [see Eq. (|18|) ]. The dimension of any 
state vector of the system is N con f. Every bosonic configuration in the respective Fock 
subspace is characterized by the set of M occupation numbers \ni, n 2 , n 3 , • • • , tim)- We can 
attribute a unique address J to each of the configurations according to the rule Eq. (1201) : 
J = J(rii,n 2 , ■ • • , njvf-i, Um)- The proposed enumeration schemes Eq. (JOJ and Eq. (120|) 
are equivalent and connected via Eq. (1T9|) to each other. They can be equally applied to 
enumerate bosonic and fermionic configurations. In other words, every many-body basis 
function written in the "fermionic style" - holes position representation - can be translated 
via Eq. (Tl9l) to the "bosonic style"- occupation numbers representation \ni,n 2 ,--- ,um)- 
However, the explicit use of the "bosonic" representation for fermions is neither appealing 
nor transparent. Therefore, once Eq. (1201) has been derived we preserve for the sake of clarity 
the holes enumeration scheme for fermions and occupation numbers enumeration scheme for 
bosons. 
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Now, we can represent the bosonic state vector in a form similar to Eq. (J7j): 

TV t N i 

1 " conj 1 * conj 

|*)= ^2Cj\J( ni ,n 2 ,--- ,n M ))= ^C^n)), (21) 
j=i j=i 

with the only difference that every bosonic configuration \ni,n2,ri3, ■ ■ ■ ,tim) has a unique 
address J which is computed according to Eq. (1201) . i.e., using bosonic occupation numbers 
representation. The index J runs over all N con f configurations. 

C. Applying operators to 

Any bosonic operator can be represented as a sum of products of bosonic creation and 
annihilation operators. If we now show that the action of any even combination of creation 
and destruction operators on a bosonic state vector Eq. (1211) leads, as in the case of fermions, 
only to re-addressing (re- indexing) of the configurations (coordinates), we can adopt all the 
ideas developed for fermions to bosonic systems. 

Let us first consider a specific, say the b 3 b 2 term which kills a boson in the second bosonic 
orbital and creates a boson in the third one. The result of the action of this pair of creation 
and annihilation operators on an arbitrary permanent is known and reads: 

b\b 2 \n 1 ,n 2 ,n^ ■ ■ ■ ,n M ) = y/n^y/n 3 + l\n 1 ,n 2 - l,n 3 + 1, • • • ,n M ). 

It acts on every many-body configuration having non-zero occupation n 2 . We can interpret 
this result as follows: The operation of any even combination of bosonic creation and anni- 
hilation operators on a permanent results in the re-addressing of this permanent to another 
one multiplied with a trivial bosonic prefactor. Clearly, all the ideas proposed for fermions 
in section [IT] can be easily adopted here for bosons. 

1. The action of one- and two-body operators 

To derive the results of the actions of bosonic one- and two-body operators on a state 
vector one can follow the strategy used for fermions in subsection III C 1[ We provide the 
final results and briefly discuss their meanings. The action of one-body bosonic operator 
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b k b q on an incoming results in |\l/ fc9 ): 

N 

1 ' conj 

\* k «) = &Kl*>=E C ?l J ( n )>> (22) 
J=l 

Cj q = Cjk q y/n^^n q + 1. 

Let us explain this expression. According to Eq. fl20|) every configuration 
|ni, • • • , rifc, • • • , n 9 , . . . , njvf) of the resulting vector |\|/ fe ^ has a unique address (index) 
J = J(ni, ■ ■ ■ , rife, ■ • ■ , rig, ... , n M ). The respective coefficient Cj q is obtained as a prod- 
uct of the J hq -th component of the incoming vector Cjk q scaled by the bosonic ^/nk^n q + 1 
prefactor. On the r.h.s. of Eq. ( 122]) we "exchange" the k-th and g-th indices due to a 
change of variables of the summation index J applied, similarly as it was done in Eq. (JHJ). 
The index J kq is computed according to Eq. (I20I) and corresponds to the configuration 
| Tii, • • • , n k — 1, • • ■ , n q + 1, . . . , n M ). This simple and straightforward methodology is ide- 
ally suitable for programming. 

The action of a general (k ^ s I ^ q) two-body b k b\bih q term on the incoming state 
vector is defined as: 

N t 

1 Y conj 

\* ksl «) = b{blb l b g \*)=J2 C J Slq \ J ( n »> ( 23 ) 

j=i 

C k / q = C Jk sig y/nk\/n~ s \/ rii + 1 yjn q + 1 , 

where, after the change of variables of the summation index J, the addresses of the incoming 
jksiq _ • • • , rife — 1, • • • , ^ s — 1? " * * j + 1) " ' " ? n g + 1, • • • , um) and resulting J = 

J(ni, • ■ ■ , n k , ■ ■ ■ , n s , ■ ■ ■ , rii, ■ ■ ■ , n g , ■ ■ ■ ,um) configurations are computed using Eq. (|2"U|) . 



2. The action of the Hamiltonian 

Using the same strategy as for the fermionic case we group together the actions of all 
one-body and all two-body operators to get the total action of the Hamiltonian on an initial 
bosonic state vector [Eq. (j2T])]. 

The one-body contributions read: 

N i 

Ay conj 

h |*> = [b\b q |¥>] \^) = °J l J ( n )) > ( 24 ) 

k,q k,q J=l 

Cj = h kq Cj q . 

k,q 
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The two-body terms are: 



N t 

J * con/ 



2 

k,s,q,L k,s,q,l J=l 

Combining the contributions from the one- and two-body terms we get the desired action 
of the Hamiltonian on a state vector |^): 

N i 

1 ' conj 

H\V) = h\V) + W\V) = C J l J ( n )>> ( 26 ) 



j=i 

c? = cj + cf, 

which concludes our constructions for bosons. 
3. Other quantities of interest 

We have seen above that the action of any bosonic operator on the state vector results in 
a new (resulting) state vector. Similarly to fermionic systems, the expectation value of the 
respective operator is immediately available as a dot-product of the incoming and resulting 
state vectors. The matrix elements of the reduced one- and two- body density matrices 
as well as the expectation value of the Hamiltonian can be obtained in a very similar way 
as have already been done for fermions in Eqs. (I13II14P and Eq. (1151) . respectively. Other 
expectation values are also amenable to this formulation. 

IV. THE CASE OF MULTI-COMPONENT SYSTEMS, BINARY MIXTURES 

In this section we generalize the new ideas for effective and efficient operations in Fock 
space proposed above for systems of spinless bosons and spin-polarized fermions to multi- 
component systems. More specifically, we show that the mapping and enumeration scheme 
discussed above can be naturally extended to more general systems. 

Let us considered binary mixtures with iV = Na + Nb particles. The mixture consists 
of Na identical particles (bosons or fermions) of type A and Nb identical particles (bosons 
or fermions) of type B. In what follows, we denote whenever needed the quantities in the 
mixture by A, B or AB superscripts. 
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The many-body Hamiltonian of the mixture has three kinds of terms: 

= £ W^A-Ak, (27) 

k,h',q,q' 

The first two terms of H^ AB "> are the A and B single-species Hamiltonians and can be read 
directly from Eq. ([T]), formally replacing in the first case the b q and b\ operators to a q and 
a k , respectively. The third term of H^ AB ^ is the interaction between the two species. We 
call H( A \H( B ) intra- and W^ AB ^ inter-species parts. 

The many-body wavefunction ^>( AB "> is a linear combination of all possible products of 
per mut at ional- symmetry- adapted configurations : 

\^ {AB) ) = i") x i™> = H c ™ ™> • ( 28 ) 

The configurations {|n}}, are either Slater determinants ([3]) or permanents fTTTj) . 

depending on whether we deal with Bose-Bose, Fermi-Fermi, or Bose-Fermi mixtures. 

In this work we make an assumption. The total number of particles of each kind Na 
and Nb is conserved. In other words, there is no conversion between the particles, i.e., the 
particles of A kind cannot become of B kind and vise versa. The Fock subspace of the A 
subsystem is spanned by all possible permutations of Na particles over M A orbitals {N A onf 
configurations) and the configurational subspace of the B subsystem by permutations of 
Nb particles over Mb orbitals {N B n j configurations). Practically, in this case we work in 
the complete configurational subspaces, i.e., the summation in Eq. ( |28l) over n runs from 1 
till N A n f, and over rh from 1 till N B n ^. More strictly, the total configurational space is a 
tensor product of two complete subspaces and any state vector of such a binary system has 
N conf N conf components. 

In the two previous sections [III and II III we have seen that enumeration schemes of the 
complete bosonic and fermionic configurational subspaces can be easily done by counting 
either holes Eq. ([6]) for fermions or particles Eq. ( l20i) for bosons. Therefore, if J a and Jb label 
configurations in the A and B subspaces respectively, then the two-component vector index 
J = (J A , J B ) enumerates all possible configurations in the total product Fock subspace. 
Now, we can represent the state vector of the mixture Eq. ( l28l) in a form where enumeration 
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of the configurations in the A and B subspaces are explicitly specified: 



N A , N B , 

conf conf 



N A . N B , 

conf conf 



\* {AB) ) = EE c 'aj b \j* c 'a'b f ) • ( 29 ) 

J A =XJ B =l J A =1J B =1 

The key goal of this section is to demonstrate that the action of the total Hamiltonian 
Eq. (12"7|) on the state vector of the mixture can also be represented and computed without 
construction of the respective Hamiltonian matrix. 



A. General fermionic system 

We consider the mixture of fermions of two different kinds A and B. Despite such an 
unusual abbreviature it is a standard system of spin-up and spin-down fermions where the A 
subsystem is associated with spin-up and B with spin-down fermions. Here, we consider the 
general, so called unrestricted case where particles are not constrained to occupied the same 
spatial orbitals, i.e., the one-particle functions of the spin-up fermions can differ from the 
respective orbitals of the spin-down fermions. Such a treatment is allowed if the Hamiltonian 
does not have terms leading to spin-flip phenomena, taking place for example in external 
magnetic fields. The Hamiltonian Eq. (|2"7j) does not have such terms. 

Let us first show that the action of typical one- and two-body terms contributing to 
the H( A \HW as well as to the parts on a basic configuration = | Jj4,(i), JbO-')) 

translates or re-addresses it to another configuration, multiplied by the respective fermionic 
prefactors. Indeed, the configuration is specified when the M£ = Ma — Na positions 
of the holes characterizing the A subsystem, i.e., the (ii, «2> • * • j*m a ) se ^i an d the = 
Mb — Nb positions of the holes characterizing the B subsystem, i.e., the i' 2 ,--- , iU) set, 
are given. Then the Ja(\) and Js(i') numbers available via Eq. (jSJ) provide unique indices 
of the two-component address of this configuration in the state vector |\I/( j4 ' B )). 

If the initial holes of the A subsystem are located at positions (h, - ■ • ,k,--- , ) i the 
action of the a\a q operator kills a particle of A kind in the g-th orbital and creates the 
particle of A kind in the fc-th orbital - we obtain the holes at - ■ ■ > ^m a )- Hence, 

according to Eq. (J6]) we get its new index J k A q of the two-component address Jb), where 
the standard fermionic prefactor has to be added to account for the correct permutation 
symmetry: 

aia q \J A (i),J B (i')) = (-l) d % J k A q ,J B ),kei,q?L (30) 
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Clearly, the actions of the two-body terms would also lead to re- addressings with some other 
fermionic prefactors: 



a\a\aia q \ JaQ), Jb^')) = j a 



J A q , Jb 



(31) 



k q 



J A lg , Jb) ,k,s G i, l,q£L 



The action of creation and annihilation operators from the A subspace changes the address 
of J A component only and does not touch the Jb ones. Similarly, the actions of the one- 
and two-body B terms of the Hamiltonian read: 



blL\j A (i),j B (i')) 



,k'q' 



J A ,J k B «'),k>ei', q '?i', 



j jk>q> 



(32) 
(33) 



,k' q ' a k , a i 

.l) d J B _1 J k B q 



j jk's'l'q' 



, k', s' G i', I', q' & i'. 



Now, we show that the action of the inter-species terms from the W^ AB ^ part also results in a 
re-addressing of the initial configuration to another one with a different fermionic prefactor: 



\JA(i),j B (i')) = (-i)^(-i)^ 



k q 



\ d .j 



k'q' 



jkq jk'q' 
J A ' J B 



k e i, q £ i, k' G i', q G" i'- (34) 



We have shown that the action of each term of the Hamiltonian Eq. ([27)1 on a general 
configuration re-addresses it to another one within the same Fock subspace. This allows us 
to conclude that the action of the total Hamiltonian Eq. (1271) on a state vector Eq. (|29|) can 
be obtained directly without construction of the respective Hamiltonian matrix: 



N A , N B , 

conj conj 



Ja=1 Jb=1 

r H(AB) _ fj(A) fj(B) r w{AB) 

^JaJb - ^JaJb ~+~ ^JaJb ~+~ ^JaJb ' 



(35) 



where C H{A) and C H(A) can be read from Eqs. (IHll9lllO|lllfT2l . and the C^j^ can be easily 
derived using Eq. 



N A , N B , 

conj conj 



= y: w^fA-Ak \* {AB) ) = E E c lT muB^m 

k,k',q,q' Ja = 1 Jb=1 

n W(AB) _ \- tjAAB) n kk'qq' 
^JaJb - 2-^1 W kk'qq' Ly J 



k,k' ,q,q' 



j JaJb ' 



jkq jk q 

C kqjk'A-lf'^-lf'B ; k ? i,q G i, k> £ i',q' G i' 

J A J B 

0; otherwise. 



f-tkk ' qq' 
^JaJb 
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Each element of the resulting vector has a unique address (Ja(i), Jb{^)) characterized by the 
two sets of holes ■ ■ ■ , q, ■ ■ ■ ,ij\,i A ) an d ■ ■ ■ ,</,•■• , z'm-b)' i- e -5 J a and Jb are obtained 
by using Eq. The value of the element C^j B) is computed as a sum of the C kq y q > com- 

A B 

ponents of the incoming state vector scaled by the respective integrals W^^ q , and fermionic 
prefactors. The address (J^ 9 , J% q ) of each of these components is obtained for every given 
set k, k,' q, q' by applying Eq. (J5j) to , i m a) and (z'i, • • • ,k',--- , i' M s). Finally, 

the expectation value of the Hamiltonian as well as of any other operator are available as a 
dot-product of the incoming and respective resulting vectors in a very similar manner as it 
has been done in the single- component case, Eq. ( floT) . 



B. Mixture of bosons 



Here, we consider the mixture of bosons of two different kinds, A and B. In section UTTl we 
have seen that operations with single-component bosonic systems can be done without rep- 
resenting the respective operators in the matrix form. Now, we demonstrate the usefulness 
of this theory and ideas to bosonic mixtures. 



The general configuration Jj = | ^(n), Jg(n')) is specified by the bosonic occupation 
numbers corresponding to the first (ni,n 2 , • • • ,um a ) and second {n'^n^ ■ ■ ■ ,n' MB ) bosonic 
subsystems. Let us first show that the actions of typical intra-species terms from the A 
subsystem lead to a re-addressings: 



a\a q | Ja(ii), Jb{^)) = Vn k + ly/n~ q 



(37) 



a\a\aia q | J^(n), Js(n')) 



a'ai 



jkq j 

J A , Jb 



jkslq j 
■J A i J B 



(38) 



Similar expressions can be obtained for the action of the intra-species terms associated with 
the B subsystem. Clearly, the inter-bosonic term acting on a general configuration translates 
it to another one, weighted by the respective bosonic prefactor: 



jkq jk'q' 
J A i B 



(39) 



Having verified that the actions of intra- and inter-species terms of the bosonic operators 
on a general configuration result in its re-addressing with some scaling prefactors, we can 
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explicitly express the result of the action of the Hamiltonian on a state vector Eq. (1291) . i.e., 
to sum up the contributions from all the terms of Eq. (J2j 



N A t N , 

conf conf 



H<4*) \*(AB }) =J2J2 \J A (n), J B (n')) (40) 

J A =1 Jb=1 

^JaJb - u JaJb ~+~ ^JaJb ~+~ u JaJb ' 

where C H(A) and C H{A) can be deduced from Eqs. (I22II23II24|25II26I) . and the C^j^ can be 
derived using Eq. fl39l) in a very similar way as done in Eq. fl36|) . 

We conclude that the operation of the Hamiltonian representing a binary mixture of 
bosons on a state vector can be performed directly without constructing the respective 
Hamiltonian matrix. 

C. Mixture of bosons and fermions 

In the two preceding subsections we have seen that effective and efficient operations 
and manipulations with state vectors of a binary mixture of bosons or a binary mixture 
of fermions are possible without resorting to the matrix representation of the respective 
operators. Actually, the final equations ( |35l) and ( )40|) for H^ AB > |v|)( AB ))) are almost identical, 
with the only exception of the enumeration scheme used to address fermionic and bosonic 
configurations. In the Fermi-Fermi case we use the positions of the holes i, i' of the A and 
B subsystems to specify the fermionic configurations and Eq. ([6]) to compute the indices of 
the two-component address of the mixture while in the Bose-Bose case we utilize bosonic 
occupation numbers n, n' to specify the bosonic configurations and Eq. (1201) to compute the 
components of the respective address. The goal of this subsection is to demonstrate, for 
completeness, the validity and applicability of the theory to a mixed system of bosons and 
fermions. 

Let A be the fermionic subsystem and B - the bosonic one. Now, to specify the gen- 
eral configuration f"j = \ Ja^), ^s(n)) one has to provide a set of the fermionic holes 
(ii,i2, ■ ■ ■ ,iM A ) an d a set of the bosonic occupation numbers (n 1; n 2 , • • • ,n MB ). The two- 
component address of this configuration is defined by the two numbers Ja(i) and J^(n) 
computed by using Eq. ([6]) and Eq. ( |2D1 . respectively. 

Clearly, the operation of the one- and two-body fermionic terms on the general configura- 
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tion changes the positions of the fermionic holes and does not affect the bosonic occupation 
numbers. In other words, only the first (fermionic) index of the two-component address 
is changed, like in Eqs. (I30|l3ip . Similarly, the action of the pure bosonic terms on a 
general configuration re-addresses only the bosonic index, analogously to Eqs. (I37|l38p . The 
Bose- Fermi interaction terms lead to the change of both fermionic and bosonic parts of the 
two- component address: 



ala q bl,b q , \J A (i), J fl (n)) = (-1)^^7+1^ J" A q , J%* ) , k G i, q £ i. (41) 



k q 



jkq jk'q 1 



Every configuration characterized by i = (ii,i2, ••• , «m a ) an d n = ( n i, n 2,-'' i n M B ) an d 
having, therefore, the address ( J^(i), Je(n)) is translated to a configuration (J^jJ^ 9 ) 
where the k-th fermionic hole was filled and a new hole at q has been created, i.e., 
(«!,••• , q, ■ • ■ , «m a ); an d simultaneously one boson from bosonic orbital q' has been trans- 
fered to orbital k', i.e., (ni, ■ ■ • , ny + 1, • • • , ry — 1, • ■ ■ , n^ B ). 

The action of any term of the Bose-Fermi Hamiltonian on a general configuration leads 
to its translation, i.e., re-addressing with some known prefactor. This allows us to apply 
the developed theory and write down the result of the action of the Hamiltonian on a state 
vector of the system as follows: 

N A t N B t 

conj conf 

#ab> = j2i: ctr: \j a ®, j B (n)> , m 

Ja=i Jb=i 

u JaJb - U JaJb U JaJb ^JaJb ' 

where the contributions C H{A) and C H( ~ B) from the actions of the intra-species Hamiltonians 
can be deduced from Eqs. (l8f9|10Mlllfl~2l) and Eqs. ( |22l23|24ll25ll26l) . respectively. The results 
of the action of the Bose-Fermi interactions CY A j^ ] can be derived using Eq. (|4T!) in a very 
similar way as done in Eq. (1361) . 



D. More components 

In the previous subsections we have shown that the total configurational spaces of 
quantum systems made of two kinds of particles, i.e., Fermi- Fermi, Bose-Bose and Bose- 
Fermi mixtures can be labeled by two-component vector index J = (J a ,Jb)- A system 
with a larger number of components can be addressed by a multi-component vector index 
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— * 

J = (Ja, Jb, Jci • • •)• Then, all the experience collected in this paper can be expanded 
to these systems as well. The only constrain is that the total number of particles of each 
kind is conserved, i.e., there are no terms in the Hamiltonian leading to particle conversion. 
We recall that, depending on the quantum statistics of the subsystem, i.e., whether we are 
dealing with fermions or bosons, we apply different schemes to enumerate the respective 
configurations. Both derived enumeration schemes, i.e., Eq. ([6]) for fermions and Eq. (1201) 
for bosons are applicable only if Na = const, Nb = const, Nq = const, .... We just mention 
here that a generalization of the presented enumeration schemes can be adopted also to 
systems with particle conversion. This issue is out of the scope of the present study. 

The most relevant conclusion is that for multi-component Hamiltonians the action of 
each term on a state vector leads to re-addressing of the configurations with some quantum- 
statistics-dependent, but simple and known prefactors. We can find the results of actions of 
each of these terms on the state vector independently. The total action of the Hamiltonian is 
obtained by summing over all the resulting vectors. Having at hand the initial and resulting 
state vectors we can readily compute the expectation value of the Hamiltonian as a simple 
dot-product of these two vectors. 

V. REMARKS ON IMPLEMENTATION 

Here, to discuss a strategy for a practical implementation of the ideas presented in all 
previous sections we refer to the single- component bosonic or fermionic system. To define 
the system means to specify the number of particles N , their quantum statistics, and the 
number of the orbitals M. The integrals hk q and Wksqi specify the Hamiltonian or any 
other operator of interest. Depending on the type of the quantum statistics the N and M 
define the length of the state vector N con f according to Eq. PJ or Eq. ffT8]) . respectively. 
It is equal to the number of the elements in the one-dimensional vector-array of complex 
numbers {Cj}j=^ representing the state vector of the system. 

Now we explain basic computational steps needed to get the result of the action of the 
Hamiltonian on a state vector. As a first step, for a given incoming state vector 
i.e., a vector-array of expansion coefficients {Cj}j^ f , the action of every pair or quartet of 
creation and annihilation operators is evaluated using Eqs. (!H|9|) for fermionic or Eqs. f)22|23p 
for bosonic systems. The result of the action of each such an operator on a state vector is 
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another state vector, i.e., as an outcome we get again a one- dimensional vector-array of 
some other complex numbers {Cj}^Z\ of the same length N con j. It is very important to 
stress that the action of each of these operators can be computed independently, implying 
an effective parallelization. For example, each available computational node is designated to 
a specific pair or quartet of creation and annihilation operators. Next, the incoming array 
{Cj}j=i is broadcasted to every node and the respective actions take place producing 
the resulting arrays. We apply Eqs. (I13II14P on each node to compute the corresponding 
element of the reduced one- or two-body density matrices as dot-products of the incoming 
and respective resulting vector-arrays. Then, by multiplying each and every elements of the 
resulting vector-array by the corresponding integral hk q or Wk sq i we get on each node the 
desired action of the respective Hamiltonian's term. Now to sum up the resulting vectors 
from all the nodes we can use an appropriate collective operation and get the desired total 
action of the Hamiltonian on the incoming state vector. Finally, using Eq. (fT5l) we easily 
compute the expectation value of the Hamiltonian (^>\H\^>) as a dot-product of the total 
resulting and incoming state vectors. 

This technique has several advantages. First of all, it does not require the evaluation of 
the Hamiltonian matrix elements in the given many-body basis set. Consequently, there is no 
need to construct, store and address these elements of the Hamiltonian matrix at all. Second, 
the elements of the reduced one and two-body matrices are immediately and naturally 
available. Third, this technique can be easily extended to three- and higher-body interactions 
potentials. Last but not the least, such a strategy implies very effective parallelization 
schemes, which are of high demands in modern computational physics. Finally, in the 
formulation of this method we did not specify explicitly the h^ q and Wksqi numbers, therefore, 
it is valid for general many-body Hamiltonians or for any other operators represented in 
the second quantization form. Consequently, this scheme can be successfully applied to 
standard real-space Hamiltonians as well as to discrete ones, e.g., of the Bose- Hubbard 
type. The derivations of the re-addressing scheme has been done for the complete Fock 
subspace spanned by permutation of N particles over M orbitals, but in principle, any 
selected subspace can be used. For example, an implementation of additional constraints on 
possible excitation patterns for fermions or restrictions on occupancies of the higher bosonic 
orbitals lead to considerable reduction of the respective configurational subspaces. In these 
cases the enumeration schemes derived above have to be modified accordingly. Moreover, 
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one has to pay additional attention to the re-addressing cases leading beyond the selected 
configurational subspace. This opens, on the other hand, a new vision on size-consistency 
issues js, 13] - it can be now explicitly considered and analyzed. 



VI. SUMMARY AND CONCLUSIONS 

In this paper we provide a novel, effective and general technique to construct the result 
of the action of any operator represented in the second quantized form on a many-body 
state vector. Within a standard framework one represents the corresponding operator in 
a matrix form and obtains the desired action of the operator by applying a matrix-to- 
vector multiplication. We have shown that the same result can be reached without even 
referring to the respective matrix elements. Considering configurations as coordinates of 
the many-body state vector, we first demonstrated that the action of any even combination 
of creation and annihilation operators on a configuration translates or re-addresses it to 
another configuration. In other words, we have seen that such an action is equivalent to 
permutation of coordinates of the initial state vector, weighted by some trivial prefactors. 
The total action of any operator, represented in the second quantized form on a state vector 
is a sum of the actions of all its terms. 

Next, for a complete subspace of the configurational space spanned by permutation of 
N fermions over M orbitals we present a simple and compact scheme to enumerate the 
fermionic configurations according to the given set of the hole's positions. Then, using the 
formal isomorphism between the configurational spaces formed by this fermionic system and 
the system of N bosons distributed over the corresponding number of bosonic orbitals we 
invent a simple and compact scheme to enumerate bosonic configurations according to the 
given set of the bosonic occupation numbers. Using these enumeration algorithms we directly 
construct the result of the action of any pair or quartet of the annihilation and creation 
operators on a state vector. Moreover, the respective matrix elements of the reduced one- 
and two-body density matrices are naturally available as simple dot-products of the initial 
and resulting state vectors. This allows us to combine the total action of the Hamiltonian as 
a sum of all its terms into simple and compact formulas, which can be directly implemented 
and permit straightforward parallelization. The proposed ideas to operate with fermionic 
and bosonic Hamiltonians have been directly extended to binary mixtures of fermions and 
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of bosons as well as to the Bose-Fermi mixtures. We have also shown that the same ideas 
can be naturally applied to systems made of a larger number of components as well. 

Finally, since the H \^f) is the basic building block appearing in the computations of 
statical properties as well as of the evolution dynamics of many-body systems, we expect 
that the implementation of the theory into the respective computational approaches will 
increase their efficiency and enable applications to systems made of larger number of particles 
than currently possible. Indeed, the proposed ideas have recently been implemented for 



bosonic systems |lOl . Illl ] within the multiconfigurational time-dependent Hartree for bosons 
(MCTDHB) 12 L and applications to multi-boson long-time dynamics in double-well [h]] 



and triple-well 



ll| traps have already been performed successfully. 
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